Metastability of a granular surface in a spinning bucket. 
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The surface shape of a spinning bucket of granular material is studied using a continuum model 
of surface flow developed by Bouchaud et al. and Mehta et al. An experimentally observed central 
subcritical region is reproduced by the model. The subcritical region occurs when a metastable sur- 
face becomes unstable via a nonlinear instability mechanism. The nonlinear instability mechanism 
destabilizes the surface in large systems while a linear instability mechanism is relevant for smaller 
systems. The range of angles in which the granular surface is metastable vanishes with increasing 
system size. 
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I. INTRODUCTION 



Granular materials display very complex dynamics due 
to their both fluid-like and solid-like characteristics EH . 
Even simple experiments can produce unexpected re- 
sults. Examples are recent experiments in which mea- 
surements were taken of the granular surface of a bucket 
of sand spun about its cylindrical axis . At low rota- 
tion rates a central region was observed in which the slope 
is significantly less than the critical slope [Q-gJ. This cen- 
tral subcritical region was conjectured to be due to the 
inertia of the flowing grains (J] since it could not be ex- 
plained in terms of a simple theory which assumed the 
surface was everywhere critical 0. 

Many models have been developed to describe the evo- 
lution of a granular surface p|-p^[. These models often 
assume that the flow is restricted to a thin layer of grains 
near the surface. Such models have been used extensively 
to study avalanches of a granular material in a horizontal 
rotating drum Jl6[ ] . However, the spinning bucket exper- 
iment is a better testing ground for these ideas since, at 
least at low rotation rates, the fundamental assumption 
of a thin flowing layer is more likely to be true jL7) . 

The purpose of this paper is two-fold. The first is 
to apply the flow models to the spinning bucket exper- 
iment. The second is to explore the source and limits 
of the metastable behavior observed in the model. In 
particular, a continuum description of the surface flow 
developed by Bouchaud et al. |l(]] and Mehta et al. |ll]] 
is used. This model includes the effect of the inertia of 
the flowing grains and is known to exhibit metastability 
via a linear instability mechanism i.e., the surface 
remains metastable for a range of angles beyond the min- 
imum angle of repose. This difference in the minimum 
and maximum stable angles was argued to be the physi- 
cal source of the Bagnold angle JToJ . 

Here I show that this model qualitatively reproduces 
the central subcritical region. I determine the mecha- 
nism creating the subcritical region and show that it is 
closely intertwined with the metastability of the model. 
The metastable surface becomes unstable via a nonlin- 



ear instability, i.e., an instability to very small but non- 
infinitcsimal perturbations. This nonlinear instability 
determines the dynamics for "large" systems while the 
linear instability is dominant for "small" systems. The 
range of metastable slopes and, hence the Bagnold angle, 
depends on the system size and vanishes in the limit of 
large systems. 



II. SUMMARY OF EXPERIMENT 

Figure |l| shows the setup for the spinning bucket ex- 
periment [Bpl] . A bucket of radius R is partially filled 
with common building sand. A conical pile is prepared at 
the angle of repose, Of — 34 ± 1°, by slowly dropping the 
plastic spheres onto the center of the stationary bucket. 
( Here we label the angle of repose 9f since, in an ideal 
Coulomb material, the surface is everywhere critical and 
the angle of repose is the same as the angle of internal 
friction The rotation rate u> is then slowly increased 
until it reaches the desired value and the resulting surface 
shape is measured. 

For an ideal Coulomb material the granular surface 
shape is found by assuming that the friction force is sat- 
urated and relating the net force on a surface element to 
the centripetal acceleration ■ The critical slope S- 

is a function of radial distance r and to, 



S-(r) 



dh 
dr 



(r/R )Q 2 - tan 6>/ 
l + (r/R )n 2 tan Bf 



(1) 



where h(r) is the local height, £1 = \J io 2 R j g is the di- 
mensionless rotation rate, g is the gravitational accel- 
eration and R Q is the bucket radius. The — superscript 
indicates that the friction force is inwards since the grains 
flow outwards. Note that S~ (r) is negative for small £1 so 
that a unstable surface corresponds to dh/dr < S~(r). 

Experimentally, the surface agrees well with Eq. pi for 
intermediate rotation rates 2 < < 5 ||,|| but less well 
for larger and smaller Q. This paper will focus on the low 
rotation rates, £1 < 2. Figure || shows the experimental 
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FIG. 1. A cut-away view of the the spinning bucket ex- 
periment. A conical pile is prepared at the angle of repose 6f. 
The cylinder is then spun about its vertical axis at rotation 
rate u> and the resulting surface shape is measured. 



surface shape obtained by Baxter for SI = 0.578 along 
with the prediction for a ideal Coulomb material using 
Of — 34° [Q. The corresponding slopes are shown in the 
inset. There is reasonable agreement between experiment 
and theory at the outer edge (r > R Q /2). However, in the 
center, r < i? /2, the slope is much less than the critical 
slope. This was thought to be because the critical theory 
neglects the effects of grain inertia B . 



III. THE MODEL 

To analyze the spinning bucket experiment we use a 
model of the granular surface evolution developed by 
Bouchaud et al. |l(| and Mehta et al. 11 1. The model 
assumes that all grains are stationary except for a thin 
layer of flowing grains on the top. We define h(r, t) as the 
local height of the immobile pile and p(r, t) as the thick- 
ness of the rolling layer. Here r is the two-dimensional 
vector giving the projection of the local position onto the 
x — y plane. 

The dynamics of the rolling layer p(r, t) is given by 



dp 
Of 



V ■ (vp) 



D V 2 p + T ({h},{p}), 



(2) 



where v = v i is the local velocity of the rolling layer, 
D Q is the diffusion constant and r o is the rate of conver- 
sion from immobile to mobile grains and V and V 2 are 
the two-dimensional gradient and laplacian respectively. 
Assuming that the grains flow outward, the conversion 
rate is: 



r„(W,M) 



( |^ - 5 (r) ) +K..V-/' 



dh 
Or 



S-(r) 



where TL is the Heaviside function. The first term in 
r o accounts for the conversion of static grains to rolling 




FIG. 2. The surface shape at Q. = ^tu 2 Ro/g = 0.578 
as given in Ref. Q. The center of the bucket is at r = 
and the edge at r/R a = 1. The dashed line is the result 
for a ideal Coulomb material with Of = 34° and the points 
are the experimental data. The inset shows the experimental 
and theoretical slopes. The experiment and theory agrees 
near the outer edge of the bucket but the slope at the center 
is significantly less negative than that of an ideal Coulomb 
material. 



grains when the slope is steeper than the critical slope 
(dh/dr < S~) and vice versa for a shallower slope. The 
pn V 2 h term causes valleys to fill and peaks to smooth. 
The last term in r o accounts for the jarring loose of static 
grains by the rotation of the bucket even in the absence 
of a flowing layer. Finally, assuming that bulk rearrange- 
ment of material can be neglected, conservation of total 
material requires 



dh 
dt 



— = -T ({h},{p}). 



(4) 



For a — the model is the same as that introduced by 
Bouchaud et al. Mehta et al. allowed for a nonzero 
a Q and also included an additional bulk arrangement 
term 

To compare with experiment we require a rough esti- 
mate of the model parameters: 

• The speed of the rolling grains v a can be estimated 
from the speed of a grain falling through a grain 
diameter d : 



I'n 



• A narrow bump of rolling grains starting at the 
center will reach the edge of the bucket in time 
t = R /v - Diffusion causes this bump to spread 
to a width Ax — V ' D t = \J D R /v when it 
reaches the edge. Since Aa; ~ y/R~o the width of 
the bump is much less than the system size for large 
R . Therefore the "large" system limit corresponds 
to 



(Ax) 2 

~rT 



Do 
v R 



< 1. 



(5) 
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The "small" limit corresponds to the width of the 
bump being of the same order as the system or 
D /(v R) > 1 @. 

• 70 AS* do/vo is the probability that a rolling grain 
jars loose a static grain as it rolls over it. A typ- 
ical excess slope of AS = 0.1 and, assuming a 10 
% probability that the static grain is converted to 
rolling, gives 7 « v /d . 

• The ratio k /^/ defines a length scale on which the 
hole filling/peak smoothing mechanism dominates. 
This is important only at very small length scales 
on the order of a few grain diameters so k w <i 7 . 

• a Q depends on the noise in the apparatus. The 
generation of rolling material by noise is assumed 
to be much less effective than the convection so that 
a Q < v - 



Equations (|2}|4|) can be rewritten in dimcnsionless form 
by measuring lengths in terms of the bucket radius R a 
and time in terms of t Q = R /v , 



dp 

dt 



-£ = -v • (f P ) + DV 2 P + r({h},{p}), 



(6) 



where 
and 



if -s-(r)H *r-7, 



dh 
Or 



-S-(r) 



«<f-S-(r) 



^ = -T({h},{p}). 



(7) 



(8) 



The dimensionless parameters are D = D /(v R), 7 = 
joR/vo, k = k /v and a = a /v . 

Assuming a large system, order of magnitude estimates 
of the dimcnsionless parameters can be obtained from the 
earlier estimates: 



D= — <1, 7= ~ -r, 

v R v d 

K d Q 0i o 

k = — « — - 7, a = — < 1. 

v R v a 



(9) 



The speed of the rolling layer is unity under this rescaling. 
As shown in the next section, the behavior of the model 
is most sensitive to 7 and less sensitive to the exact value 
of D and n as long as D <C 1 and rc/7 <C 1 . 

An especially interesting feature of this model is that 
it displays metastability and hysteresis. Bouchaud et al. 
performed a linear stability analysis by balancing the rate 
at which a small bump of the rolling layer is convected 
downhill with the rate at which the bump is amplifies 
and diffuses |lfj. They found that the bump affects the 



behavior uphill only if the slope exceeds the critical slope 
by an amount larger than 



AS = S~ - — > 



dh 
Or 



Dole 



1 

D7' 



(10) 



The metastable behavior was interpreted as the physical 
source of the Bagnold angle, i.e., the excess angle be- 
yond the angle of repose at which a static sandpile first 
becomes unstable. This dynamical interpretation of the 
Bagnold angle is very different from the usual mechanical 
interpretation of this angle JhJ . 



IV. NUMERICAL RESULTS 
A. Spinning Bucket 

To analyze the spinning bucket experiment, Eq. 
were numerically integrated in polar coordinates with az- 
imuthal symmetry. Using the experimental parameters, 
6 f = 34°, d a « 0.25 mm and R a = 12.5 cm j§, we 
obtain 7 = R /d = 500 and rc/7 = d /R = 1/200. As- 
suming a "large" system and small vibrations, we choose 
D = 0.001 < 1 and a = 10~ 8 . A third-order Adam- 
Bashford method was used ]19| with mesh-size Sx = 0.002 
and time-step St = 5 x 10 . Zero slope boundary condi- 
tions are applied at r = and zero curvature boundary 
conditions at r — 1. Any rolling material reaching the 
edge is converted to static. 

An initial conical pile is prepared at the critical angle 
and the rotation rate f2 is increased instantaneously from 
il = to £1 — 0.6. Since the critical slope decreases with 
f2 the initial surface is unstable and, as shown in Fig. 
||, an avalanche is generated. At t « 0.125 (Fig. ||a) a 
buildup of flowing material is visible at the edge of the 
bucket. This buildup generates a steep local slope in the 
static pile which causes the static material uphill to fail. 
The failure then propagates uphill reaching the center of 
the bucket at t 0.2 (Fig. ||b). The static pile beneath 
the flowing layer approximately follows the critical curve. 
The flowing layer then convects to the edge of bucket 
(Fig. pfc) leaving behind a central subcritical region in 
which in which the slope is much less than the critical 
slope (Fig. 

The integration was repeated with different values of 
D, 7 and k to test the robustness of the results. The 
central subcritical region is qualitatively the same for all 
7 > 400 and D < 0.005 becoming more pronounced for 
larger 7 and smaller D. For smaller 7 and larger D (7 < 
200 and D > 0.02) the final surface follows closely the 
critical surface. Varying the ratio k/j or a did not have 
noticeable effect as long as they were small but nonzero. 

A comparison of model results in Fig. |^ with the ex- 
perimental data in Fig. || shows that the essential fea- 
tures of the experiment are reproduced. In particular, 
there is a subcritical region in the center and a region 
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FIG. 3. The evolution of a conical pile at f2 = 0.6. The 
dashed line is the critical surface and the shaded area is the 
flowing layer. In (a) a flowing region has formed at the edge 
of the bucket. This failure zone propagates uphill reaching 
the center in (b). The static surface under the flowing layer 
closely follows the critical surface. The flowing layer convects 
downhill (c) and the final surface (d) has a central region with 
a slope less than the critical slope. 



near the edge which follows the critical surface. How- 
ever, a more detailed analysis shows two important dif- 
ferences. First, the flowing layer is much larger than the 
few grains assumed in the model. Second, in the exper- 
iment, the rotation rate is slowly ramped up from zero 
rather than changed instantaneously. Hence the numer- 
ical results corresponds to the actual experiment only if 
the final surface is independent of the initial state, or, if 
the initial conical surface is metastable up to an fl close 
to 0.6. 

To mimic the experiment more closely we repeat the 
integration while slowly ramping the rotation rate from 
n = to f2 = 0.6. We monitor the excess slope be- 
fore an avalanche occurs, i.e., the amount the slope of 
the metastable surface exceed the critical slope, and the 
magnitude of the subcritical region after an avalanche, 
i.e., the amount the slope of static surface falls below the 
critical slope. 

The conical surface is metastable as f2 is increased from 
zero until the first avalanche occurs at f2 » 0.2. The dy- 
namics of this avalanche is very similar to that for the 
instantaneous jump to ft = 0.6. The main differences are 
that the excess slope before, and the magnitude of the 
static slope after the avalanche are both much smaller 
for the instantaneous jump to f2 = 0.6. As is in- 
creased further the critical slope decreases and becomes 
shallower than the new static slope. Eventually a sec- 
ond avalanche occurs at O as 0.35. The process repeats, 
leading to a third avalanche at ss 0.45 and a fourth 
at Q, s» 0.55. After each avalanche a subcritical region 
is observed with a much smaller magnitude than for the 
instantaneous jump. This indicates that the magnitude 
of the subcritical region after an avalanche increases with 
the excess slope before the avalanche. This was confirmed 
by instantaneously changing the rotation rate to different 
values of f2 and observing the subcritical region after the 
resulting avalanches. 



B. Two dimensional sandpile 

The numerical results indicate that the central subcrit- 
ical region is closely related to the metastable behavior 
of the model. To understand the factors determining 
the limits of metastability we consider the simpler case 
of a two-dimensional sand-pile in a stationary container 
(fi = and a = 0). This eliminates the complications of 
the cylindrical geometry and the constant production of 
flowing grains by the bucket rotation. 

We start by preparing a static critical surface. We first 
produce a surface at the critical slope: dh/dx = — tan Of 
along with a small uniform rolling layer. Here x is the 
linear position measured in units of the system size L 
with the high end of the pile at x — and the low end 
at x = 1. This constant slope surface is not stationary 
due to the finite value of k. To obtain the static critical 
surface we integrate Eq. (Iq-Bj) in one dimension until the 
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FIG. 4. The two-dimensional sandpile immediately after 
the failure zone has propagate all the way to the top. The 
static layer below the rolling layer follows closely the critical 
slope (dashed line). However, as shown in the inset, the slope 
of the static layer (solid line) is slightly less steep than the 
critical slope (dashed line). 



initial flowing layer is completely converted to static. 

Once this static critical surface is obtained, we desta- 
bilize it by tilting the surface through an excess slope 
AS. We also assume the tilt produces a uniform flowing 
layer Ap. To determine the characteristics of the result- 
ing avalanche, we used a large AS = 0.1 with the same 
parameters as we used for the spinning bucket 7 = 500, 
D = 0.001, K/7 = 1/500 and Ap = 10" 8 . The mesh-size 
was Sx = 0.001 and time-step 5t = 2 x 10 -7 . The result- 
ing avalanche has the same features as the avalanche in 
the spinning bucket: 

1. The avalanche is induced by the growth of the 
rolling layer at the foot of the pile. 

2. This buildup of the rolling layer creates a valley in 
the static pile. The large local slope of the static 
pile causes the static material uphill to fail. 

3. The failure zone propagates uphill until it reaches 
the top of the pile. As shown in Fig. ||, the sur- 
face of the static pile below the rolling layer follows 
closely the critical surface. However, a more de- 
tailed examination (inset of Fig. ^J) shows that the 
surface of the static pile is slightly less steep than 
the critical slope. 

4. Since the propagation of the failure zone is faster 
than the convection of grains down the hill, the 
total height of the pile (static + flowing) is the same 
as the initial static pile. Therefore, the built-up 
flowing layer at t* , the time the failure zone reaches 
the top of the hill, is 

p(x, t*) « h(x, 0) - h(x, t*) = AS(l-x). 



This is confirmed in Fig. ^ which shows a trian- 
gular shaped flowing layer with p(0, t*) = AS and 
p(l,t*=0. 

5. The rolling layer flows downhill leaving a subcritical 
region at x < 0.5. The surface for x 3> 0.5 follows 
closely the critical surface except at the very edge 
where there is a buildup of the static pile Q . 

Since the magnitude of the subcritical region after an 
avalanche depends on the excess slope before, we need to 
explore the limits of metastability. The linear stability 
analysis predicts that the surface is stable until the ex- 
cess slope exceeds AS* = l/(Dj) 10 1. However, several 
characteristics of the avalanche indicates that a differ- 
ent mechanism may destabilize the surface. First, the 
instability always begins at the bottom of the pile, in- 
dicating that the system-size is important. Second, the 
avalanche is generated when there is an increase in the 
local steepness of the static pile. The linear stability anal- 
ysis neglects the effect of the flowing layer on the static 
pile. 

To check if the linear stability dominates the dynamics, 
we determined AS**, the minimum excess slope required 
to initialize an avalanche, as a function of the different 
parameters. Fig. ^|a shows AS** as a function of the am- 
plification rate 7. Although both the linear stability pre- 
diction and the observed AS* are proportional to I/7, 
the observed values are a factor of 20 smaller the linear 
prediction. 

A more serious conflict is shown in Fig. ||b where AS"* 
is plotted as a function of D. The numerical results 
clearly deviate from the linear stability prediction that 
AS** oc 1/D for small D. Hence, for large systems, i.e., 
small D, the linear stability mechanism is not relevant to 
the surface evolution in this limit. On the other hand, for 
large D, AS* ~ 1/D so the linear instability mechanism 
may be relevant for small systems. However, for large 
D, the subcritical region disappears and the final state 
follows closely the critical surface. 

Fig. ||c shows AS as a function of the rolling layer 
Ap. Here AS* oc — log(Ap) so that the metastabil- 
ity of the initial granular surface depends on the am- 
plitude of the initial perturbation. In particular, a small 
but non-infinitesimal rolling layer is required to induce 
an avalanche indicating that the avalanche is initiated 
through a nonlinear instability mechanism. 



V. THEORY AND ANALYSIS 

A. Nonlinear instability mechanism 

To understand the metastable behavior observed nu- 
merically two effects must be considered: First, the flow- 
ing layer is amplified as it flows downhill, and second, 
the amplification of the flowing layer changes the slope 
of the static pile. 
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FIG. 5. (a) The critical excess slope AS* as a function 
of I/7. The dashed line is a straight line fit showing that 
AS ~ I/7. (b) ASa* as a function of 1/D. For small D (or 
large systems) AS is only weakly dependent on D. For large 
D, AS* ~ 1/D indicating the linear instability may initialize 
the avalanches in small systems, (c) AS* as a function of the 
initial rolling layer Ap. The dashed line is a straight line fit 
showing that AS ~ — log(Ap). 



To do so we modify an argument given by Bouchaud 
et al. for the amplification of a bump of rolling grain as 
it flows 113]. Consider a critical surface tilted through 
an initial excess slope AS with an uniform flowing layer 
Ap. Neglecting any change in the static pile, Eq. (||) 
shows that the flowing layer will grow exponentially as 
it flows downhill. For small D, the flowing layer will 
approximately maintain its shape so 



p{x,t) 



if x < t, 

Ap e^ Ast ]fx>t. 



(11) 



In terms of the dimensionless variables both the radius of 
the bucket and the speed of the grains is unity. Therefore 
for t > 1 all flowing material reaches the bottom leaving 
the origin static state intact. 

This argument indicates that the granular surface is 
stable for all AS. However, this is the case only if the 
change in the static pile can be neglected. The change in 
h(x,t) for small D is 



Ah(x.t) = h(x,t) - h(x,0) 



dt' T({h},{p}) 



- f dt' 7AS p(t') 
Jo 



(12) 



Substituting Eq. ( |Tl| ) for p{t') and integrating to t > 1 
gives (since p(x, t) = if t > x), 



Ah(x,t) = -Ap (■ 



s 7ASa; 



(13) 



Therefore the flowing layer generates an additional excess 
slope of 



AS'(x) 



dAh 

dx 



jAS Ap e 7ASa; 



(14) 



If AS" is not small there will be a positive feedback mech- 
anism. The increased steepness makes the flowing layer 
grow faster, which in turn, generates an even larger ex- 
cess slope. Once this positive feedback mechanism builds 
up, the increased slope causes the material uphill to fail 
and the failure zone propagates up to the top of the hill. 
Since AS'(x) is largest at x = 1 the avalanche must start 
at the bottom of the pile in agreement with numerical 
observations. 



Using Eq. (14) the minimum initial excess angle re- 
quired for an avalanche occurs when A«S"(1) w 1 or 
7 AS* Ap e T AS * w 1. Solving for AS* gives 



AS* 



log(Ap) + log( 7 AS*) _ logAp 

7 7 



(15) 



where we assume Ap -C 1. This result is in agree- 
ment with our numerical results of the previous section: 

AS* - I/7, AS* log(Ap) and AS* approximately 

independent of D for small D. 

A more intuitive understanding is obtained by rewrit- 
ing the metastability criteria in terms of the original 



G 



unsealed parameters. Using 7 = j L /v and Ap = 
Ap /L to rewrite Eq. [IB] in terms of the dimensionful 
variables we have 



Therefore the slope of the static pile after the avalanche 



AS* 



j L "° \ L 



•log 



(16) 



Therefore, according to this model, the range of 
metastable slopes and hence, the Bagnold angle, depends 
on the system size L Q . Furthermore AS* vanishes as 
L — > 00 showing that the metastability vanishes for in- 
finite systems. This is because it takes time t Q = L /v 
for the flowing material at the top of the hill to reach 
the bottom. During this time the flowing layer is am- 
plified by a factor e 7 ° AS *°. The larger the system the 
more the flowing layer is amplified and the more likely 
an avalanche is initiated. 



B. Central subcritical region 

We now return to the origin of the central subcritical 
region observed after an avalanche. Assume an avalanche 
is initiated at the bottom of the pile. For large 7 and 
small D, we empirically find that failure zone propagates 
uphill much faster than the rate at which grains are con- 
verted downwards. We find that the static pile below 
the flowing layer is slightly less steep than the critical 
surface (see Fig. ||). Since the slope of the static pile 
is less than critical, the built-up rolling layer is continu- 
ally converted to static as the layer flows down the hill. 
The flowing layer is convected away from the top of the 
pile first so this conversion occurs longest at the bottom, 
leading to final surface that is shallower than the critical 
surface. 

This argument can be made more quantitative by 
approximating the built up rolling layer as a constant 
p(x,ti) = AS where t\ is the time the failure zone 
reaches the top of the pile. As a further approxima- 
tion, we also assume that the underlying static pile has 
a constant slope slightly less steep than critical dh/dx = 
— tan Of + Sd where Sd is the deviation of the static pile 
from the critical slope. For small D we can neglect the 
change in the shape of the flowing layer as it is convected 
downhill: 



p(x,t) 





AS e 



7 ASt 



if t • 
if * 



t\ > x, 
t\ > x. 



(17) 



The change in the static pile after the flowing layer has 
passed is 



is 



Ah{x,t) 



dt' r({h},{p}) 



Substituting Eq. ( |i~7] ) for p(t') and integrating to t — 1\ > 
1 gives the total change in h, 



dh 

Ox 



- tan 9 f 

- tan 6 f 

- tan 8 f 



dAh 



- Sd 

- Sd 

- ^Sd AS e 



dx 

jS D AS e " ySDX 

— jSox 



(20) 



where we neglect the second term since 7 is large. This 
result is in agreement with our qualitative picture. In 
particular, the deviation of the static pile from the critical 
surface is largest at small x. 



VI. SUMMARY 

We have applied a model of granular surface flow de- 
veloped by Bouchaud et al. |l(]] and Mehta et al. |ll]] 
to the spinning bucket experiments of Baxter ||Q . The 
model qualitatively reproduces the central subcritical re- 
gion observed in the experiment at low rotation rates. 
The subcritical region occurs when a metastable surface 
becomes unstable via a nonlinear instability mechanism. 
The nonlinear instability is due to the amplification of 
the rolling layer as it is convected downhill. This ampli- 
fied layer causes the static pile underneath to become 
steeper which in turn causes the flowing layer to be- 
come even larger. This positive feedback mechanism ini- 
tiates the avalanches in "large" systems. We numerically 
determined the excess slope required to destabilize the 
metastable surface as a function of system parameters 
and showed that it agrees with our analytic arguments. 
In particular, the model predicts that the excess slope, 
and hence the Bagnold angle, depends on system size and 
vanishes in the limit of large systems. Lastly we showed 
how the nonlinear instability leads to the central subcrit- 
ical region via a conversion of rolling grains to static as 
they roll downhill. 
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dt' jS D p(x,t'). (18) 



Ah(x) = AS (1 



) 



(19) 
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